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Abstract 

Without resorting to special treatments for each individual test case, the ID and 2D CE/SE 
shock-capturing schemes described in Part I [1] are used to simulate flows involving phenomena 
such as shock waves, contact discontinuities, expansion waves and their interactions. Five I D 
and six 2D problems are considered to examine the capability and robustness of these schemes. 
Despite their simple logical structures and low computational cost (for the 2D CE/SE shock- 
capturing scheme, the CPU time is about 2 /rsecs per mesh point per marching step on a Cray 
C90 machine), the numerical results, when compared with experimental data, exact solutions or 
numerical solutions by other methods, indicate that these schemes can accurately resolve shock 
and contact discontinuities consistently. 
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1. Introduction 


The al nlity to generate finely resolved, shocks and contact discontinuities at correct 
locations without introducing numerical oscillations is an important requirement of an 
accurate solver for high speed in viscid flows. Few established solvers l if indeed they exist) 
can meet this requirement consistently without resorting to special treatments (such as 
tuning some ad hoc parameters) for each invidual case. The purpose of this paper is to 
show that, despite their simple logical structures and low computational cost, the LE> and 
2D C'F/SF shock- capturing schemes described in Part T [L| generally can meet the above 
requirement consistently with no special individual treatmen ts required. 

Tt was shown in [2] that the ID OF,/ ST. shock-capturing scheme is highly accurate in 
solving Sod’s shock-tube problem. Recently, its accuracy in solving TTarten's shock-tube 
problem was also evaluated against other schemes by Ratten cl- /)/. [3|. They concluded 
that “Tl i e result (see Fig. T) is quite remarkable, considering that the internal details 
of the Riemann fan are never used and, indeed, are never needed, because the entire 
Riemann fan is contained within the region of in tegration ." In this paper, the accuracy 
of the present LD scheme will be further evaluated using additional test problems. The 
results to be presented will further substantiate another conclusion reached in [3|, i.e., L Tn 
fact, solutions produced with this scheme are- broadly comparable to conventional MUSC'Tj 
schemes " 

As the first step, three standard test problems involving an infinitely I eng shock tube 
are used to evaluate the performance of the LD CTL/SF. shock-capturing scheme. Flow phe- 
nomena associated with these test problems include shock waves, con tact discontinuities, 
expansion fans and their interactions. In the second step, two test problems involving a 
shock tube with two closed ends are used to evaluate the capability of the LD C'F/SF 
shock-capturing scheme to simulate shock reflection accurately. 

Finally, six 2D test problems are used to evaluate the performance of the 2D C'F/SF 
shock- capturing scheme. The p hen 0111 en a associated with these problems include complex 
shock-wave interactions, and multiple reflections of shock waves from solid bodies. One can 
conclude from the numerical results that the C'F/SF method maintains its advantage over 
the traditional methods in simplicity, generality and accuracy even in a multidimensional 
space. 

T 11 the following sections, the boundary conditions used will be specified prior to the 
description of the numerical results. Note that the symbols and notations used in Part T 
( L] are inherited in Part TT. 
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2. Boundary Conditions 


2.1. Boundary Conditions for ID Tost Problems 

Note that numerical simulations of the first three LD test problems to be discussed 
in Sec. % are- termin atet l before the region of nonuniform flow reaches the right and left 
boundaries of the computational domain. As a result, the boundary conditions used in 
these problems are- immaterial as long as they are compatible with the uniform nature of 
the solution at each boundary. For simplicity, the steady-state boundary conditions are 
imposed at both boundaries. The solid-wall boundary conditions imposed at the closed 
ends of the shock tube involved in the last two ID test problems are discussed here. 

Consider Fig. L. Let (i) x\ — .tv, Jax with 2 Y, being a given constant and j a half or 
whole integer, liil the mesh point with j — r> — 0 uof- belong to the set fl, (iiij it 1 ■ 0 be a 
half or whole in teger, and [ivl the right and left boundaries of the shock tube be located at 
x — 2 s T and x — 2 -__ iT , respectively. The solid-wall boundary conditions at each boundary 
will be constructed by assuming that, at any time f, the flow fields to the right and left 
of the boundary are the mirror images of each other. This assumption coupled with the 
definitions of n L , and given in Tvj. id. L l of [2| implies that, at any time I and for any 
x, (i) 


On 


On, 


: x, f ) — x . i — x, Cl and ' (j*. ; j*, h — x : : — x, Ll (2.L) 


Ox ' 


Ox 


where m — L,3; and liil 


Oil 2 

i j 2 1. j* j j, h — — i j 2 1_ je : ^ — x, i l and — — j, fl — 


Oil** . 


x, f) (2.2) 


Tjet Li) — L/2, i/2, d/2, if.it is a half integer; and iiij — 0, L,2, if jt is a 

whole integer. Then the numerical boundary conditions corresponding to Tvjs. (2.L) and 
(2.2) are 


l l 


— l n. 


i 11 . 

1 - l/2 


an ( l i 


i 11 . 
j .i 


- j* - L/2 


I ! J i 


m-L,3 \2.X\ 


and 


I u-> l 


— —I i I 


and 


r , 'i 11 

■ ~ r - - A : l/2 


— I I 


l'2.di 


respctively. Note that, according to Fig. L, the boundary points ( . jt, f>) ■/ 11 if and it 
are related by either of the relations li) and iiij given preceding Tvj. (2.3). 

Tjet j't l)e a half integer, then ( . jt,0j fl. As a result, the initial data should be 

specified at the mesh points l j, 0), j — — jt, — .it L, ,jt~ L_it- Note that the data must 

be consistent with the requirement that velocity vanishes at j — . it- 

Let it be a whole integer, then ( . jt, 0) ■/ fl. As a result, the initial data should be 
specified at the mesh points lj,0j, j — — jt L/2,— jt 3/2 jt — 3/2, jt — L/2. 
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Given the 1 above initial data and the boundary conditions Tvjs. (2.3) and (2.4), the 
marching variables at all l j, u) fl witl i .(;, ' j ' — j't and r> ' ■ 0 can be obtained by using 

the LD GTl/fjn shock-capture scheme, i.e., Frjs. (2.32) and (2. bo) of Part T fl|. 

2.2. Boundary Conditions for 2D Test Problems 

The steady-state oblique shock problem suggested by Yee and others [4| is one of 
the 2D test problems to be discussed in Sec. 3. This test problem will be used as the 
prototypical case in the following discussion. 

Tl, e computation domain and the shock locations (.if and FF\ are depicted in Fig. 2. 
Tl le lower boundary is a solid-wall. Assuming y — L.4, the eract Filler solution to this 
problem is: 

la) Tii the region . i/?f. 


n-2.Q, ft- L.O, y>- L .0/ L A (2.o) 

[ b ) Tii the region .\FFD, 

ii- 2.(5193, i* - -0.50(532, ^-L.TOOO, ji-L.o2$2 (2.(5) 

(c) Tn the region FCF, 

m - 2.40 La, ?’-0., f> — 2.b£T2, p - 2.0340 (2.T) 

Tl ie mesh used in the current numerical calculations is depicted in Fig. 3. Again a 
mesh point ( j, h 1 , u I i." fl L (i.e., u — 1/2, 3/2, o/2, . . .1 is marked by a solid circle; while a 

mesh point (j, k\ »j (i.e., u — 0, L , 2 , ) is marked by a hollow circle. The mesh 

is a special case \ b — 0) of that depicted in Figs. T— 9 in Part T [L]. Note that only the 
mesh points i." fG are- present at the inflow and outflow boundaries. Moreover, a mesh 
point and the corresponding marching variable will be identified by the time-level number 
u, and two uew mesh indices r and which are- given in Fig. I as a pair of integers 

enclosed in parentheses. Note that, for the mesh points fl L , r — L,2,3 R, R L, and 

£ — 2, 4, (5 2 S. On the other hand, for the mesh points fib, r — L,2,3, , R , R L, 

and a — L , 3 , b 2.S — L,2S L. 

With the above preliminaries, the initial and boundary conditions can now be spec- 
ified. At all mesh points l r, s, 0 1 , it is assumed that (ij m — L,2,3,4 are- evaluated 
usiii K r G- I 2.o I, and l ii l 

m-L, 2,3,4 [2.A\ 

Furthermore, for u — L, 2, 3, , the above conditions (ij and (ii) are- applied at all mesh 

points at the inflow boundary A/? (see Fig. 3). 

At the upper boundary .if 1 , for all r> — L/2, L, 3/2, 2, , ii'i are- evaluated using 

Dj. ( 2.(3), and (ii) n . and :j ' are evaluated using Tvj. (2.S). 
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The- solid-wall boundary conditions at DC will be constructed by assuming that, at 
any time I, the flow fields below and above DC are the mirror images of each other. By 
using the definitions of m — L, 2, 'l, given in Tvj. (fi.L) of Part T fL| and the fact that 
y — 0 at any point on DC, it can be shown that the last assumption implies that 


— y, ( i — T-.-, (j*, y, f }, rn — L,2,d, and ;jj( j 4 , — y, ( i — — uji j,y, f) 


l 2.0 1 


0'n . 0\i t-i . . On Tyi . On tvt . i n m ■ n i ■■■.■ 

— — I j*, — y, 1 1 — — — |j%y, h and — — I j 4 , — y, f I — = — |j%y, h, hi — L,2,d i2.L0i 

O? ' (Jj? Otj ' t.iy 


an c l 


On j . On 3 . On s . On 3 . 

-^r— I J-, — y, I j — I.J\y, f.l and — — l j 4 , — y, 1 1 — — l j 4 , y , 1 1 

d* d" ()<j Lhj 


l'2.LL I 


Consider the mesh depicted in Fig. Then the numerical boundaty r conditions corre- 
sponding to Fifjs. (2.0)— [2. L L) are 


Mr. Mrs 


(2. L 2 i 


an ( l 


(’n-nr) R j and (i^i' ) X . l 3 0^1')^ f>! L,2,d 1.2. L^.l 




respectively. According to Fig. 3, the range of « in Ivjs. (2. L2j— (2. Ld) is dependent on the 
time level h. Tjet (i) S f — S if $ is even; and and (ii) S* — S — L if S is odd. Then (i) 

— d,S, L2 2S f if h — L/2, 3/2 and (ii) „-s — L,o,0, ,2S 1 L if r> — L,2, 

Furthermore, by using Fvj . (C.d) of Part T (L| with b — 0, it can be shown that Tvjs. (2. L:Vl 
and (2.L-1) are- equivalent to 


' s ' R ■ L ,J 


(' i ^Sr^ (' i ™Sr. (j^c) ' 


(2. Lo) 


R 


and 


W 0^)r. 


I n -rs- ^ 


3<T 


; R,s 


I 2. Lfi) 


respectively. Equations (2.L2), (2. La) and (2. Lfi) are- the boundary conditions at the lower 
wall (a solid-wall). Tn other words, the marching variables associated with the mesh points 
below the solid-wall will be determined using these equations. 


At the outflow boundary CD, for any h — L,2,3 and r — L,2,3 R, we assume 


that 


/. v 1 _ /, y i L /- 

L. '* 111 ■ I K'*™) <7 ' 




■ L 


-o. 


hi — L,2,3,d 
c>i — L,2,3,d 


(2. IT) 

(2. Lsi 
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m — 1,2,3, "3 


i 2. ISd 


ail < L 


i ( 1 


T1 J i 


)" 
■ T* 


,2$ 


By using Tvj. (C.d) of Part T fL] tritli b — 0, it ran hr shown that Ftjs. (2. LSj and (2.19) 
are ecjuivalen t to 


an ■: l 


L < 


• ll-L/2 


2 




L < 


■ n-l/1 


2 0 

i -n-tj 1 

iJ TMr J „ ’ 



m — L,2,3,d 


I 2.20 1 


l 2.2 L I 


Tl iils the marching variables at the outflow boundary Trill hr determined using Ffjs. (2. LT), 
(2.20) and (2.2 L). Note that, according to thr numerical results to hr presented in Sec. 3, 
thr outflow boundary conditions arr non -reflecting in naturr. 

With thr aid of thr ahow initial and boundary conditions, thr marching variables at 
any timr Irvrl can hr drtmninrd hy using thr 2D CF/SF shock- capturing sclirmr. 

Tn thr oblique shock prohlrm drscrihrd ahow, only a horizontal solid Trail is present. 
Tn other 2D test problems to hr drscrihrd later, hotli horizontal and vertical solid-trails may 
hr present. As trill hr sliotrn immediately, imposing thr solid-trail boundary conditions 
at a vertical trail over a uniform mesh similar to that depicted in Fig. 3 is slightly mote 
complicated than that at a horizontal trail. 

Consider thr mesh depicted in Fig. dial. Let TIC and CD hr solid-trails. Note 
that, given any exterior mesh point ( R L,.<s, »} that lies imntctlialfty helotr DC, erne 
can find an interior mesh point that lies at the same time level and also is the mirror 
image (relative to DC) of the exterior mesh point. As a result, the solid-trail boundary 
conditions Fqs. (2. L2(— (2. Ld) can he imposed. Contrarily, given any exterior mesh point 
(r, 2S L, i?l, — 0, L,2,... that lies j \mmctlial-cly to the right of CD, one can iiof. find an 
interior mesh point that lies at the same time level and is the mirror im age (relative to 
CD) of the exterior mesh point. As a result, the mirror image conditions for the vertical 
trail CD, i.e.. 



,iT * 

r»i — 

L,3,d, and 


-Ml* 



and 



m ~ L , 3, d 


and 


/ \ in 

/ \ 11 

and 

/ \n , 

vn 

L ” 


L ” “( 

. 1 J — J - 1 


(2.22) 

(2.23) 


(2.2-3) 


can not he used directly as the solid-trail boundary conditions. 

Ttro approaches can he used to overcome the ahow difficulty. Tn the first approach, 
the following approximate forms of Fqs. (2.22)— (2. 2d) are- used as the solid-trail boundary 
conditions: 


-64)^ tn — L,3,d, and (W)^,. t “ “64)^ (2.2o) 
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i L ' f>! U i. 2.2(5 ,i 


an c l 


L ’ al1t ^ ( :i ^)^. L ( ;j 2i , ) t ,^J 


(2.2T1 


Here f>! — air- evaluated using the kiinmi marching variables at- the 

mesh point- (r,25’,f> — L/21 with the- aid of the first-order Taylor's expansion. 

Tn t-lii-’ second and generally mote accurate approach, the present matching procedure 
is applied over a dual mesh. It was explained in Sec. S of Part T that a dual mesh can he 
considered as the combination of tiro staggered space-time meshes, i.e., mesh L (Fig. -If all 
and mesh 2 (Fig. 4(b)). As shown in Fig. 4(b), a mesh point of mesh 2 is also identified 
by tiro spatial indices r and a, and the time-level number t>. Furthermore, a mesh point 

( r, *, t? I of mesh 2 is marked by (il a solid triangle if t> — L/2, :{/2 , o/2 and (iil a hoi low 

triangle if t> — 0, L,2 

Let (il r and a be any pair of whole integers and t> be any whole or half integer such 
that ( r, a, t>) represents a mesh point of mesh L, and liil C — t> . L/2. Tl len, according to 
Figs. -1(a) and 4( bl, ( r, a, t>* 1 represents a mesh point of mesh 2 at the »'th time level with 

its spatial location being identical to that of I r, a, i> I. Furthermore, (il for i> — 0, L,2 

(r, 2£ L,i>) is a mesh point of mesh L while l>,25, r> j is a mesh point of mesh 2, (iil lor 

— L/2, 3/2, 0/2 l>,2,5' L,ulisamesh point of mesh 2 while ( r, 2. S', i>) is a mesh point 

of mesh L, and (iiil fir — 0, L/2, L,3/2 the mesh points ( r, 2$ L,i>) and (r, 25', r> j 

lie at the same time level and, relative to CD, are mirror images of each other. As a 
result, the solid-wall boundary conditions at CD can be imposed using Fqs. ( 2.22 )-( 2.24), 

with — 0, L/2, L, 3/2 Note that the marching variables associated with meshes L 

and 2 are- now coupled through these boundary conditions. Also note that the Fqs. (2.231, 
(2.2-11, (2.2(31 and (2.2T1 can easily be converted to the versions associated with the ( C , tj 1 
coordinates by using Fq. (C.41 of Part T fL| with b — 0. 

Imposing the conditions Fqs. (2.221— (2.211 directly requires the use of a dual mesh. As 
a result, it has the disac Lvan tage of doubling computational cost. However, the extra cost 
no longer becomes an issue if the use of dual meshes is mandated by other considerations. 
As will be explained in future papers, this indeed is the ease for many practical applications 
including those requiring the use of unstructured triangular meshes. 
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3. Numerical Results 


3.1. One-Dimensional Problems 

As a preliminary, a subtlety related to the specification of numerical initial conditions 
rrill be discussed here. As an example, t re consider a shock-tube problem defined by the 
conservation law Ffj. (2. L9j in Part T [L| and the exact initial conditions: at I — 0, 


■.U-f 
Ol — 


i fj*' ■(/ 

r~, if J ■ ■ d 


(3.L) 


Here li) m — L, 2, 3, are the flow variables that appear in the LD Puler equations 
Tvj. [2. LA) of Part T fL], and (ii) and F~ (F^, -f- F~j, r>i — L, 2, 3, and d are given 
con stall ts. 

To proceed, consider a computational mesh (see Fig. L j with the properties that (ij it is 

a whole number, and (ii) a mesh point (j‘,0) fl if and only if j i." {. L/2, . i/2 Ijt — 

L/2)}. For each ij’,0) i." fl, the line segment joining the two space-time points (j* : . aj"/ 2,0) 
is part of SF( j, 0). As a result, the space- time fluxes passing through the above line segment 
can be evaluated using either the exact or numerical initial conditions. The resulting two 
values are identical if and only if 












IJ 


(j*,0}dj*, 


m — L,2,3 


(3.2) 


Because of its flux-based nature, accuracy of the C'F/STI Method generally will sutler 
[ patLlmtatty i if the exaet initial eomlitlonfi are not eoitfuiiious) if the numerical initial 
conditions specified does not satisfy Tvj. (3.2). 

Tjet d — JVj, with jj C {0, . L, .2, , . ijt, — L)}. Then (jj, 0) is the mid- point of the 

line segment joining the mesh points ( jj . L/2, 0) fl. As a result, Tvj. (3.2) is satisfied if 
the numerical initial conditions are-: (ij 


an d I ii I 


f 

if J — JJ 

L/2, jj 

3/2, . . . 

ti )';! — ;j TV1 i> J :,0) — e 




: J 1 

if J — J J 

- L/2, jj 

- 3/2, . . . 

■ fl fyi . . . . 

- o 

-■ 0 * 

for j — 

- L/2, . 3 

•■'2 


I 3.3 1 


i 3.1 1 


Note that (ii with (n^ jS being defined by Tvj. (3.3), Tvj. (3.2) is satisfied regardless what 
the values are assigned to (ii^. )S, and (ii) the impact on numerical results as a result of 
replacing Tvj . (3.1) with other choices for l generally fades away quickly after a few 

time steps. 
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T/t with j'j { . L/2, . 3/2, , . (j't L/2)}. Then (jj, 0} represents a mesh 

point- i." fl. Furthermore, according to Flq. (3.L), and 0\i ^ /V?r are not defined at the 1 
location of this mesh point. As a result, the initial conditions Flqs. (3.3) and (3.d) are no 
longer applicable. However, it can be shown easily that FIq. (3.2) is satisfied if, for any j 
with lj',0) fi, we have (i) 


T-1 l i fi.b ^ I _J JlF 

r-)/2, if j jj 


(3.5) 


and (ii) 

l>,^-)° — 0 ( 3.(3) 

Note that fi) with (ij^ jS being defined by F5q. (3.5), FLq. (3.2) is satisfied regardless what 
the values are- assigned to and (ii) the impact on numerical results as a result of 

replacing F5q. (3.(5) with other choices lor ( , )' : : ' generally fades away quickly after a lew 
time steps. 

Newt we will com pare the total n timber of the space- fime mesh points involved in a 
CFI/SFI simulation and that involved in a simulation using a typical regular- mesh scheme 
which requires only one marching step to advance by a time period Ad Fet (i) T. be the 
spatial dimension of the computational domain, (ii) T be the total simulation time, and 
(iii) and f>o be the whole integers with I. — faAr and T — »oAb Then, for the staggered 
mesh un der eon side-ration (<o — 2_Jt in this case), the total number of the space- time mesh 
points involved in a 0FI/SF1 simulation (excluding those that lie on the initial line and 
spatial boundaries) is (2<o — L) ■ f>o. On the other hand, the corresponding number for a 
regular-mesh scheme is ( (V, — l) ■ nv,, i.e., about one half of that involved in the 0F1/SF1 
simulation. The difference is due to the fact that it requires twomarching steps to at Ivan ce 
over each time period A I in the CF/ftF method. 

Furthermore, we shall assume in this section that (ii p, u , jj, e, d/j, and fi — L.Tl denote 
the mass density, velocity, static pressure, sonic speetl, shock wave Mach number and 
specific heat ratio, respectively, anti (ii) C*F7. denotes the maximum value of ( |u | ejAl/AJ* 
reach et l t luring a numerical simulation. Also, for each test problem, (i) a mesh similar to 
that depicted in Fig. L will be used in numerical simulation, and (ii) the exact solution (if 
it exists) will be plotted (as a solitl line) along with the numerical solution in every figure. 

Also, note that o (see Flq. (2.(55) in Part T [L|j is the only adjustable parameter in the 
ID 0T1/SF1 shock- capturing scheme. Generally, numerical wiggles near a discontinuity are 
satisfactorily suppressed by setting rt — L. However, in some cases, complete suppression 
may require the use of a larger value of o (e.g., o — 2,3). 

With the above preliminaries, the numerical results of five LE> test problems will 
be presented. The first three involve shock tubes of infinite length. They are- selected 
to demonstrate the IE) 0F5/SF, shock-capturing scheme's robust capability lor capturing 
shock and contact discontinuities as well as lor simulating expansion waves. The last two 
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test prol ile ms in volw shock tubes wit-1 1 closed ends. They are- selected to dem on strata that 
the- present solver is an accurate tool to study the end effects of shock- tube flows. 

3.1.1. Sjogrccn Problem 

Tl iis problem is taken from [5| . The initial conditions at i — 0 are: 


{ [ L.0,0.1, —2.0), i ■ ■ 0.5 

[3.T) 

( L. 0,0.1, 2.0), if j*- -0.5 

The CTi/Sn solution at ( — 0. L [Fig. 5) is generated assuming o — L, — 0, j! 4 iT — L, 
to ~ 100 [i.e., aj 4 - 0.0 L), and Ai 1 - 0.002 iCFT. — 0.55). 

Th e initial velocity discontinuity causes two rarefaction waves to propagate in opposite 
directions, leaving in between a region of high vacuum. It was mentioned in [(51 that several 
Godunov-type schemes failed in this problem due to extremely low pressure in the middle 
region. Tt can be seen that the present solution agrees very well with the exact solution, 
without showing negative pressure in the middle region . Tts accuracy is comparable to 
that obtained by Xu ef. at. [(51 using a gas-kinetic scheme with (q — 200. 

3.1.2. Slm-Osher Problem 

Th eiiiteraction ofa moving shock [.!/., — 3} with a sinusoidal density wave is examined 
in this problem [T] . The initial conditions at I — 0 are: 

I" [3.*5T, 10.333,2.(529), if * • ■ -4 
< [3.*) 

j [ L 0.2 sin 5 j», L. 0,0.0), if j* ‘ ■ —4 


The OFi/SF, solution at I — L.S [Fig. (5 ) is generated assuming o — L, j 4 -^ — —5, j iT — 5, 
f 0 - *00 [i.e., aj* - 0.0 L25) and A i - 0.00 15 iCFT. - 0.5*2). 

Tl iis problem does not haw a known exact solution. Several upwind schemes haw- 
been used to solve this problem to com pare- their ability in resolving the peaks appearing 
in the solution [*]. The present solution is coin parable to those obtained in [*] by using 
the TVD L and TVD2 schemes with fo — *00. 

3.1.3. Merging of Two Shock Waves 

The phenomenon considered in this problem involves two shock waves propagating 
to the right in an infinitely long tube with a shock of — 3 behind a weaker shock of 
1/j — L.5 [pp. Li L- Lid of [91). The initial conditions at I — 0 are-: 


f [ T. L*23, 25.10 1(5, 3.*12(55), if * ■ ■ -2.0 


} [ L.*(j2L,2.l5*' 
[ L.0, L. 0,0.0), 


— { [ L.*(52L, 2.15*3, 0.*2L(5), if —2.0 ■ ' j 4 ■ ' L.5 [3.9) 

if L.5 ■ ' j 4 


NASA/TM— 1998-208844 


10 



A OE/SF, simulation is carried out- assuming o — 2, ■n , —j h — —3, r ; y — T, ("o — LOO (i.e., 
aj' — O.L) avkL aI — 0.0 L35 [OFT. — 0.*2). The numerical t^^nl At- i — 0.f5T5, L.L205 and 
Lfi2, shown ill Figs. T— 0, respectively, agree very ’roll with the exact solutions. The two 
shocks mil alii separated until l — L. L205. After tliAt, one shock catches up with Another 
And they merge in to a stronger shock. The lAtter propAgAtes to the right At a higher speed, 
leAviiig behind it a contAct surface And a left-moving expansion fan (see Fig- &). 

Note that there is a sin aI I pressure overshoot near the shock in Fig. 0. As shown in 
Fig. L0, practically identical numerical results at i — L.(32 are obtained if o — 2 is replaced 
by rt — 3. However, the overshoot is no longer present in Fig. L0. 

3.1.4. Woodward-Colella Problem 


Tl iis problem, concerning the interaction of two blast waves in a tube with closed 
ends, was proposed by Woodward and Colei la f L 0] . It has no known exact solution. The 
initial conditions at I — 0 are: 


1(1.0,1000,0), if r ■ 0.L 


ip,p, u ) 


< ( L .0,0.01,0), if O.L ■ ' r • ' 0.0 

[ I L.0, LOO, 0), if 0.0 ■ ' r 


(3.L0) 


The two ends are, respectively, at r — 0 and r — L where the reflecting boundary conditions 
are- imposed. 

The C'E/SF. solution at l — 0.03* (Fig. LL) is generated assuming o — L, (\-j — *00 
(i.e., aj* — 0.00 L2.o ) and A I — L.25 x L0 - " 1 [CFT. — 0.3524). The flow field at i — 0.03* 
(containing three contact surfaces and two shock waves) shows that the contact surfaces 
are more smeared than the shock discontinuities. The current numerical results are- at 
least as accurate as those f L l| generated by using the AUSM ' , Roe, Van Leer, AUSMDV 
and AVSM ' -w splitting schemes (see Fig. 12). 

3.1.5. Waves in a Shock Tube with Closed Ends 


The flow pi lenomena studied by this problem (see Fig. L3) include (ij the reflection 
of a shock wave and an expansion wave, and (iij the interaction of the reflected waves and 
contact surface. The tube has closed ends at £• — 0 and jp — L, respectively. At I — 0, a 
diaphragm located at jp — 0.25 separating two gases at different conditions is burst l pp. 
205-20* of [01 ). With the initial conditions 


\p,p, a ) — 


(20.0,20.0,0.0), if jp 
( L.0, L.0, 0.0), i fj* 


0.25 


0.25 


(3.11) 


the breaking of the diaphragm creates a shock wave and an expansion wave separated by 
a coil tact surface. 
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To demonstrate the present solwr's capability to maintain its accuracy after a long 
running time, the present simulation with o — 2, Ar — 0.0 L and A i — 0.00b is carried out 
for a period that is long enough to shorn - successive reflections and interactions of shock 
waves, expansion craves and contact surfaces. When I — 0.00, the waves haw not yet 
reached the ends as shown in Fig. M. When t — O.b, the left-moving expansion waws 
have been reflected from the left end while the shock waw is still moving to the right (see 
Fig. Lob. At ( — O.d, the shock wave has already completed its reflection from the right 
end and the reflected shock wave is moving to the left (see Fig. LfJ). This reflected shock 
interacts with the right-moving contact surface to create a left-moving transmitted shock, 
a right-moving contact surface, and a tight-moving shock waw that, when reaching the 
right end, is reflected again. At I — O.oSo, there exist in the flow field two shock waws, one 
contact surface, and a region of expansion waws (see Fig. LT). Numerical results plotted in 
Fig. Ld-LT correctly display the flow phenomena described in Fig. Lb at the four design ated 
time instan ts. 

3.2. Two-Dimensional Problems 

Til this part, six problems in volving shock reflection, diffraction, and interaction are 
solved using the 2D CTI/STL shock- capturing scheme. The first is a steady-state problem 
while the rest are- time-dependent. A steady-state solution is obtained as the converged 
solution of the time-marching procedure- [L2]. Numerical results are compared with the 
exact solutions, experimental data or numerical solutions obtained by using other methods. 
Tn all of the following numerical simulations, o — 2 is assumed throughout the entire- 
computational domain. Because o is the only adjustable parameter in the present scheme 
( see Sec. fi in Part T f L] ), the .■same .scheme js applied a t a// id ferior me.sh poid f„s regardless of 
the presence of bow disco;; f.jdujf.je.s. Also the approximate reflecting boundary conditions, 
which were discussed in Sec. 2 and can be implemented without using a dual mesh, will 
be imposed on a verf-ieal solid-wall boundary unless specified otherwise. 

To paw the way for the following presentation, a further discussion of the space- time 

mesh depicted in Fig. b is in order. According to Fig. b, at each time level u — 0, L,2, , 

there are S L staggered columns of mesh points (marked by hollow circles) with each 
column containing R L mesh points. Thus there are- \ S L ) X ( R L) mesh points at 
each of these time levels. Furthermore, because two neighboring columns are- separated 
by a distance ir while two neighboring mesh points in any column are- separated by a 
distance 2h, we haw S’ — ITVir and R — Tf/i 2h l where IF and 7 f are- the width and I light 
of the computational domain, respectively. 7f the drsf, fJ:e f-li/rd, the fifth , — eoJumd.s 
were moved upward a disf-adee h , then the mesh poidfs marked by hottom’ circles form a 
regular Carfesiad .spa ibal mesb wj(J: S' add R mesl; In fervals m the r- add y- djreefiod.s, 
respectively. As a result, the mesh formed by the mesh points marked by hollow circles 
will be referred to as a 5 X mesh. 

Similarly, at each time level u — L/2, b/2, , there are- S’ staggered columns of mesh 

points (marked by solid circles) with each column containing R L mesh points. Again two 
neighboring columns are separated by a distance ir and two neighboring mesh points in 
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any column arc- separated by a distance 2 h. In this paper, the mesh formed by these mesh 
points mil also be referred to as a S X R mesh. Zii gonotal, wtrat-fdSo.oo of licw its ootumno of 
titcfili points aw positioned, a spa (fa! mesh eoveniig a wetanerutat eomput-ational domain of 
width IF and height Tf v/HS be referred to ax a f 1 T^/'n ■ j x l Tf /[ 2h I j tnesh if two neighJ?e>t-ing 
eotumns of mcol: pom to aw .separated by a distsnee ir and two neig}d noting me.sb poin lo- 
in Aiix- eoSumn aw separated by a fdistanee 2 ft. 

Fora dual mesh referred to earlier, there are fro set of mesh points at one time level. 
Tf each set forms a (TTVir) X I Tf /[ 2h 1 1 mesh, then the dual mesh mil be referred to as a 
dual idr, dr) x I TT/\ 2h ) I mesh. 

From the above analysis and the fact that it requires fro marching steps to at Ivan ce 
by a time period At in the CF/SF method, the total number of s/iarr-Z-jme mesh points 
involved in a 2D CF/SF simulation is approximately equal to \2Tf At l x S x R, i.e., about 
twice that of a 2D single-step regular-mesh simulation if each simulation uses a X R 
mesh and both have the same values of Ai 1 and total simulation time T. Note that, in the 
special case that a dual S X R mesh is used, the total number of mesh points involved in 
a CF/SF simulation is approximately equal to i/lT/Af I X S' X R. 

3.2.1. Oblique Shock Problem 

The computational domain, mesh structure and initial /boundary conditions used in 
the current simulations of this problem were described in Sec. 2. A numerical simulation 
is carried out using a (30 X 20 mesh with At — 0.0 L. The resulting steady-state pressure 
contours and the pressure coefficient O p (— 2lj)/j) (>5 — L l/l yd/”.. I with — 2.0 and 

p.-,-. — L.0/L.-1 being the in flow Mach number and pressure, respectively) at >j — O.o are- 
plotted in Fig. IS, where the solid line represents the exact solution. The im prow-men t 
on shock resolution by using a finer L20 X dO mesh can be seen in Fig. 19. No numerical 
oscillations are detected near either the incident or the reflected shocks, and the computed 
Op agrees very well with the exact solution. 

3.2.2. 2D Supersonic Plow Past a Step 

Consider the supersonic channel flow of db, — 3.0 past a step depicted in Fig. 20, 
which is a standard benchmark problem in the literature. The flow exhibits complicated 
phenomena such as Mach stem, slip surface, expansion fan, and their interactions and 
reflections. Tt was used to test Flarten’s TVD VLTLC scheme [13], Oiann akouros and 
Kamiadakis's spectral elemen t-FCT ( flux- corrected transport) method f Ld], and Van Leer's 
ultimate conservative difference scheme f La ]. Tt was also used by Woodward and Colella f L 0] 
to compare the accuracy of different numerical methods in handling a shock discontinuity. 

Note that the upper comer of the step is the center of a rarefaction fan and hence 
is a singular point of the flow. According to Woodward and Colella [10], unless special 
numerical treatments are- applied near the comer of the step, the computed solutions would 
be seriously affected by large numerical errors generated just in the neighborhood of this 
singular point. Specifically, these errors cause a boundary layer to form just above the step 
in the wind tunnel. A shock then interacts with this boundary layer, and the qualitative 
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nature of the flow in the- tunnel is altered more or less dramatically, depending upon 
the difference scheme and the mesh used. 7t- wit! .he tmmecitat-ety that- satJkfaet-or}' 

tiumetleat zotutiotiG rad be obtauied by the pteoen t .scheme without employing fipeelat 
treatments at the upper eotnet of the step. 

Tl le mesh used in the current simulation is also depicted in Fig. 20. Note that no 
mesh point is placed at the singular point at the upper comer of the step. The initial 
conditions are- set to he the free stream conditions. Furthermore, the constant free stream 
conditions are- imposed at the inlet while the nonreflecting boundary conditions Fcjs. I 2. IT), 
(2.20) and (2.2 L) are- imposed at the exit. In addition, the reflecting (solid-wall) boundary 
conditions are- imposed at all other boundaries. 

To show the improvement in flow solutions with decreased mesh spacing, the density 
con tours of the solutions obtained by the present solver with (30 X 20, L 20 X 40 and 2-10 X SO 
meshes are- shown respectively in Fig. 2L. Note that the values of A t used in the above 
computations are- identical to those used in f L 0] , i.e., 0.00 To, O.OOo and 0.002a, repectively 
[OFT. — 0.S). From Fig. 2L, it is seen that the Mach stem, triple point, slip surface, 
expansion fan at the comer, and the interaction between the reflected shock and the 
rarefaction waves are- accurately simulated in the present solutions. Note that an alternate 
simulation in which the dual-mesh reflecting boundary conditions Fcjs. ( 2.22 j— ( 2.2-1) are 
imposed at the vertical step wall yields almost iden deal results. 

3.2.3. Blast Wave Problem 

A blast flow field generated at an open-ended cylindrical shock tube, which was de- 
scribed in the experiment by Schmidt and Duffy [ 1(3), has been studied by Wang et at. f IT] 
using an axisymmetric version of the present solver. The cylindrical shock tube con figu- 
ration and the initial conditions are depicted in Fig. 22. An unsteady flow field is created 
at t — 0 by the sudden removal of a diaphragm at the lip of the tube, which separates a 
compressed fluid in region 2 inside the tube from the surrounding stagnant fluid in region 
L The direct contact of the high and low pressure regions results in a rarefaction wave- 
propagating back into the tube and, in the meantime, a shock wave blasting away from 
the tube lip into the ambient open space. The coin parison between experimental results 
and OF/SF solutions presented in [IT] shows that the unsteady blast wave development 
was captured very well in the numerical solutions. 

Consider the planar-flow version of the above problem where it is assumed that the 
configuration of the planar-flow version is such that its cross section is identical to that of 
the axisymmetric- flow version. In other words, the tube in the original problem is replaced 
by two slabs in its planar-flow version. A numerical simulation of the new version was 
carried out using the 2D CF/SF shock- capturing scheme. The computational domain is 
depicted in Fig. 23 with the understanding that the axis of symmetry is represented by 
the lower boundary. For illustrative purposes, some CF/SF mesh points are- also shown 
schematically in the same figure. 

The computation was performed over a L(3D X L 20 mesh with A i — 0.002 a. The non- 
reflecting boundary conditions are imposed at the inlet and outlet, while the reflecting 
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boundary condition); arc- imposed at the- lower and upper boundaries, and at all slab ’rails. 

To s 1 1 ow the- time history’ of dor: development, numerical solutions at eight time levels 
are shown in Figs. 24 and 2a, in which the pressure and density contours ranging from 

0 to o.S3 with a constant interval of 0.040 are- plotted. The sequential plots reveal that 
as the blast wave initiated from the open end of the region between tiro slabs propagates 
to the right, a vortex is developed at the lip of the slab real I , and it mows downstream 
with an ascending motion. When the blast waw reaches the upper wall, it is reflected as 
shown in the plot at ( — L.O msec. In the meantime on the axis of symmetry a normal 
shock is formed ahead of the vortex and mows slowly in the downstream direction, while 
a jet shear layer is created at the lip of the region between two slabs. At I — L.o msec, the 
portion of the blast waw that is reflected from the upper wall is shown to mow toward 
the vortex. After passing the vortex, the blast waw becomes curved and keeps moving 
forward to interact with the normal shock below, as shown in the plots at I — L.T msec. At 

1 — L.O msec, the flow pattern reveals that as a result of the interaction, the blast waw is 
broken into two parts while several new vortices are created. More complex flow patterns 
are- shown at ( — 2.L and 2.3 msec, describing further reflection and interaction of shock 
waves and vortices. 

Despite the difference between planar and axisymmetric con figurations, the computed 
flew fields agree extremely well with those shown in the shadowgraph pictures of the test 
results [Lfl] taken at ( — 0. 199(5 and 0.493T msec, respectively. 

3.2.4. Diffraction of Shock Wave down a Step 

Two experimental cases, described in fp. 14$, LA], about the diffraction of a plane 
shock waw down a step, are simulated here. The computational domain (0 r, y ' 3.o) is 
depicted in Fig. 20 with the upper boundary’ being a solid wall. The shock waw is located 
at the comer of the step at I — 0. The reflecting boundary’ conditions are imposed at the 
upper and lower walls as well as all step walls. The non-reflecting boundary conditions 
are- imposed at the inlet and outlet. Tn the first case in which a weak shock wave with 
U,~ L.3 is considered, the numerical simulation is carried out using a L40 x TO mesh with 
A I — 0.0 L. The computed solutions at ( — 0.42, L.L2, and LA2 are- also shown in Fig. 20, 
in which the density’ contours are plotted from 0.S to Loo with an interval of 0.02o. A 
comparison of the computed solutions with the first set of the photographs shown on p. L4S 
of f Is] (reproduced here as Fig. 2T) reveals that there is general agreement in flow pattern, 
except for those phenomen a induced by’ the viscous effect. The observed shock waws, slip 
lines and vortices are- very’ well captured in the numerical solutions. 

Tn the second case in which a stronger shock wave with — 2.4 is considered, the 
numerical simulation is carried out using a LSD X TO mesh (0 j* 4.0 and 0 y 3.o) 
with A I — 0.0 Do. The computed solutions at I — 0.2 To, O.STo, and L.3To are- shown in 
Fig. 2S, in which 30 density’ contours are- plotted with the contour values increasing from 
the minimum of O.o with a constant interval of O.OOo. The flew patterns exihibited in 
Fig. 2S are more complicated than those in Fig. 20. Tn large measure they’ are- also similar 
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to those observed in tin.’ second sot of the photographs shown on p.LdA of [ LAj (reproduced 
hero- as Fig. 29). 

3.2.5. Shock Rcflcclicm from a Dust Layer 

Flore a practical problem of shock reflection from a dust layer is studied. Following 
the wedge model described in f L 9) , we consider a plane shock moving to the right with 
Mach number IG, — L A L toward a wedge whose surface is inclined at angle as shown 
in Fig. 30. Square protuberances of size T. f 2 are placed at equal distances 7. apart on the 
surface to simulate dust particles. The common origin of the two coordinate systems (j,y) 
and (y,i/) is situated at the tip of the wedge, with the jd and t/ axes being parallel and 
normal to the wedge surface, respectively. 

As depicted in Fig. 3 L, the computational domain l —0.0 T.O and 0 '<J / 7. 

d) contains seven protuberances. The front of the incident shock thus makes an angle 
with the </ axis. At ( — 0, the computational domain is divided into two flow regions by the 
shock front that intersects the jd-axis at Jf'/T. — — O.d. Standard stationary atmospheric 
conditions are- assumed in the region to the right of the shock front, while constant fluid 
conditions with dJ, — L A L are assumed in the other region. Reflecting boundary conditions 
are- imposed at all solid walls, while non -reflecting boundary conditions are implemented 
at both the inlet and outlet, and on the part of jd axis with —0.6 ' jd/£ ' 0, through 

which waves can move freely. On the upper boundary ( </ / 7. — d), where the reflected 
waves have not reached before the end of all simulations, numerical values are assigned 
ahead of and behind the plane shock according to the exact solution. Computations were 
carried out for — 20°’, 30°, and dO 0 , respectively, using a 300 x 1^0 mesh with Af — 0.0 L, 
where I is made dimensionless using 7. as the reference length and the speed of sound in 
the undisturbed region as the reference speed. 

Til order to show a clear com parison between the experimental and computed results, 
the CTL/STl solutions are plotted in the ( r, yj coordinates through a coordinate transforma- 
tion . First, to show the unsteady evolution of wave patterns resulting from the reflection of 
shock waves over the dusty wedge, computed density contours at four difTerent time levels 
for — 30° are plotted in Fig. 32. The effect of varying wedge angle on wave pattern 
can be observed in the density contour plots shown in Fig. 33 lor — 20° at I — 3.A 
and in Fig. 3d for — d0° at i — 3.0, when the incident shock wave is standing at the 
upper right corner of the dth protuberance. The sell Keren photographs taken from f L 9) are 
reproduced in Figs. 36— 3T to show representative wave patterns for the cases <?„. —20°, 30° 
and d0° at different instants. Tn these photographs, (i) Model Fi and Model F) represent 
the laboratory models with 7. — A and 2 mm, respectively, and (iij T L denotes the triple 
point generated by the reflection of shock waves from the First protuberance. The location 
of T l in the c direction is indicated by the numerical value of jf/T. in each figure. Tt is 
seen that, as the incident shock wave moves forward, a compression wave is reflected from 
each protuberance and an expansion wave is generated from its back. Gradually, the- indi- 
vidual compression waves accumulate to form an envelope (* r and a stronger compression 
wave C 1 J, (see Fig. 3d). For the cases with — 20° and f?„, — 30°, the developments of 
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waw patterns are almost the same, while for <?,,. — d0°, a. kink point A" appears as shown 
in Fig. 3T. A com parison between experimental and numerical results indicates that the 
photographed wave patterns are- eorreetly eaptured in the OF/SF solutions. The dose 
resemblance between Fig. 33 and Fig. 3 a and that between Fig. 32(d) and Fig. 3(3 in terms 
of both wave and vortex structures are- dearly recognizable. 

3.2.8. Implosion/Explosion of Polygonal Shock Waves in a Box 

The problem concerning the implosion /explosion of a polygonal shock wave in a square 
box studied in [20] is investigated here. Not only the early stage of the implosion /ex plosion 
process, but also its later development, which was not studied in [20], are simulated here. 
Tt furtl ler demonstrates the robustness of the current Fu ler scheme in its ability to model 
multiple shock reflections and interactions. 

All simulations are- carried out using lil a dual 2d0 X 2d0 mesh covering a square box 
(—2 j-, t/ ' 2 l , and [iij a ( *FJ. number — 0.0. The reflecting boundary conditions are 

imposed at the four sides of the square box with the understanding that the more accurate 
dual-mesh reflecting boundary conditions Fqs. I 2.22 l— l 2. 2d) are used at the vertical walls. 

The initial shock wave configuration is a regular polygon. Tt is assumed that (i) the 
polygon shares with the square box the same geometric center (located at (0,0)), (iij one of 
the vertices of the polygon is located at (0, 0.$v’3)’ and (iiij there is a low pressure region 
inside the polygon with a pressure- ratio of L0 across the shock. Xote that, as a result of 
(ij and (iij, the vertices of the polygon are points on the circum lere-n ce- of the circle that 
has a radius — 0.£v"3 suid is centered at (0,0). 

As the first step, the early flow field is studied for three cases in which the initial 
shock wave configurations are an equilateral triangle, a square and a regular pentagon, 
respectively. The computed pressure and density contour plots at different time levels 
are shown in Figs. 3$ and 30. Wave patterns similar to those shown in Figs. L-o of 
[20], obtained using a TVD method on a 3o0 x 3o9 mesh, are- clearly observed in the 
C'F/SF solutions, displaying detailed features such as Mach stems and polygon-shaped 
flow discontinuities. 

As the second step, the implosion /explosion of a hexagonal shock wave is simulated 
until the second implosion of the shock wave is observed in the box. More complex flew 
phenomena can be seen in the pressure and density con tour plots of Figs. 40 and 4 L , 
including the reflections of shock waves, shock-shock interaction, and shock-contact surface 
interaction. Tt is interesting to note that the shape of the contact surface centered at the 
center of the box remains unchanged even after the passage of shock waves. Further 
development of the implosion /explosion process can also be simulated using the same- 
solver. 
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4. Conclusions and Discussions 


By coni paring the current numerical results with experimental data, exact solutions 
or numerical solutions generated by other methods, it has been shown that the LD and 2D 
CF/SF shock- capturing schemes can accurately resolve shock and contact discontinuities 
consistently. Furthermore, it has also been shown that the present schemes are genuinely 
robust, i.e., uaSskc tnatiy other fihoek-rayitui-ma; .sebesne.s, fJiei \t aorutary are achieved with- 
out resorting to .specia/ frea tticen is; for eaeh UuUvUiust ca.se. 

Because of their simple logical structures and totally explicit nature, the present 
schemes are- also highly computational efficient. As an example, consider a vectorised code 
implementing the single-mesh 2D CF/SF shock- capturing scheme. For a 300 X L 20 mesh, 
the C'PU time on a Cray COO required to execute LSfi marching steps \T — LSfi X l A t /2 l ) 
is only L4 seconds, i.e., about 2. Lb /isecs per mesh point per marching step. 

A disadvantage of the present schemes is that the total number of space- time mesh 
points involved in each single-mesh [dual-mesh l CF/SF simulation is about twice [four 
times) that of a typical single-step regular-mesh simulation if each simulation uses a S x R 
mesh and both have the same values of S, R, At and total simulation time T. However, 
this disadvantage can be coin pen sated for by other advantages the present scheme has, 
such as higher accuracy and lower computational cost per mesh point per marching step. 
As an example, consider the test problem discussed in Sec. 3.2.0. The TVD results given 
in [20] are generated using a 359 X 359 mesh while the CF/SF results are generated using 
a dual 240 x 240 mesh. Assuming that the same CFT. number is used in both simulations, 
then the value of a/ used in the CF/SF simulation is 350/240 times that used in the 
TVD simulation. Tt can be shown easily that the total number of space- time mesh points 
involved in the CF/SF simulation is only about 20% more than that used in the TVD 
simulation. The slight ( lisa: Ivan tage of the present scheme can be further coin pen sated for 
not only by its possible lower computation cost per mesh point per marching step, but 
more importantly, by its ability to simulate the implosion /explosion process long past the 
early stage simulated by the TVD method. 
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Figure 1: A staggered space-time mesh with the spatial boundaries atj = ±jb (jb = 2). 



Figure 2: The computation domain and the shock 
locations of a steady-state oblique shock problem. 
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Figure 3: The spatial locations and the mesh indices (r, s) of mesh points used 
a steady-state oblique shock problem (R = S = 4), 
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Figure 12: Upwind solutions of the Woodward-Colella problem. 
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Figure 13: Waves in a shock-tube with closed ends. 
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Figure 14: The CE/SE solution and the exact solution of the waves 
in a shock-tube with closed ends problem ( t = 0.09). 
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Figure 15: The CE/SE solution and the exact solution of the waves in a 
shock-tube with closed ends problem (t = 0.3). 
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Figure 16: The CE/SE solution and the exact solution of the waves in a 
shock-tube with closed ends problem {t = 0.4). 
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Figure 17: The CE/SE solution and the exact solution of the waves in a 
shock- tube with closed ends problem (t = 0.585). 
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Figure 18: Pressure contours and pressure coefficient at y = 0.5 of the 
oblique shock problem (60x20 mesh). 
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Figure 19: Pressure contours and pressure coefficient at y = 0.5 of the 
oblique shock problem (120x40 mesh). 
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Figure 20: Geometry and Grid distribution of the 2D supersonic flow past a 
step problem. 




Figure 2 1 : Density contours of the 2D supersonic flow past a step problem 
generated using 60x20, 120x40, and 240x80 meshes. 
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Figure 22: The initial conditions and geometry (cross section) of a cylindrical shock tube 
for the blast wave problem. 
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Figure 23: The computational domain and mesh-point distribution of the blast wave 
problem (planar-flow version). 
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Figure 24: Pressure contours of the blast wave problem at eight different time levels. 
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Figure 25: Density contours of the blast wave problem at eight different time levels. 
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Figure 26: The computational domain and density contours at three different 
time levels of the diffraction of shock wave down a step problem (first case). 
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Figure 27 : Experimental results of the diffraction of shock wave down a 
step problem (first case — three different time levels). 
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Figure 28: The computational domain and density contours at three different time 
levels of the diffraction of shock wave down a step problem (second case). 
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Figure 29: Experimental results of the diffraction of shock wave down a 
step problem (second case — three different time levels). 
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Figure 30: Shock moving past a wedge with a dust layer. 
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Figure 31: The computational domain of the dust layer problem. 
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Figure 32: Density contours for the dust layer problem (6 W = 30°) at four 
different time levels, (a) 1 = 0.5, (b) t = 1.75, (c) t = 3, (d) t = 4. 
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Figure 33: Density contours at t = 3.8 for the dust layer problem ( 0 W = 20°). 
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Figure 35: A schlieren photography for ( 0 W = 20°). 



Figure 36: A schlieren photography for ( 0 W = 30°). 



Figure 37: A schlieren photography for ( 0 W = 40°). 
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t=0.39 t=0.484 t=0.694 


Figure 38: Pressure contours for implosion/explosion in a square box with 
different initial shock wave configurations. 

(a) an equilateral triangle, (b) a square. 

(c) a regular pentagon. 
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Figure 39: Density contours for implosion/explosion in a square box with 
different initial shock wave configurations. 

(a) an equilateral triangle, (b) a square. 

(c) a regular pentagon. 
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Figure 40: Pressure contours for implosion/explosion of a hexagonal shock in a 
square box. 
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Figure 41: Density contours for implosion/explosion of a hexagonal shock in a 
square box. 
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